Bose statistics and classical fields 
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Classical fields counterpart of the ideal Bose gas statistics in a trap is investigated by performing 
calculations in the canonical ensemble. There exists the optimal cut-off which allows to match the 
full probability distribution of the condensate population by its classical counterpart. Universal 
scaling of that cut-off with temperature and dimensionality is derived. 



At the end of XIX century equilibrium thermodynam- 
ics of the electromagnetic field, the theory of a so called 
black body radiation, suffered for thepersistent ultra- 
violet divergence. It was Max Planck [1| who cured the 
problem introducing a concept of light quanta later called 
photons. With the help of photons the energy carried by 
high frequency degrees of freedom was tamed. 

More than 20 years later Einstein [3], stimulated by 
preliminary ideas of Bose 01 introduced statistical prop- 
erties of the ideal gas of indistinguishable massive par- 
ticles that are now called bosons. In particular it was 
shown that entirely due to the indistinguishability, when 
cooled, bosons tend to assemble at the ground state of 
the binding potential at relatively high temperature. 

This phenomena called Bose-Einstein condensation 
was finally realized in the seminal experiments in 1995 
Of course realistic theory had to depart from the ideal gas 
and was to account for interactions. While a well known 
Gross-Pitaevskii equation combined with small excita- 
tions - the BogoHubov approximation |5i| - forms a good 
basis of analysis near absolute zero temperature, many 
experiments require a method valid for higher tempera- 
tures. Several groups proposed a classical field approxi- 
mation as an effective approach to nonzero temperature 
Bose gas 0|. This is like traveling the road from Jaynes 
to Max Planck backward in time -ignoring the granular 
character of atoms and stipulating that amount of atoms 
may be changing in a continuous way. Understandably 
equipartition of the energy between the modes leads to 
the ultraviolet divergence just as it was for classical radi- 
ation. Hence, the authors introduce a short wave length 
cut-off to cure the problem. Different authors and some 
of them differently in different papers choose the value of 
the cut-off 0,0]. 

It is the purpose of this Rapid Communication to show, 
based on exactly soluble models of the ideal gas, that 
there exists the optimal cut-off that allows to match the 
full probability distribution of the number of condensed 
atoms by its classical fields counterpart. We give a scaHng 
of this cut-off with temperature and dimensionality of the 
binding potential. Of particular interest is the population 
of the last state retained in the optimal classical fields 
approximation. It depends on dimensionality but not on 
temperature. 

We start our analysis with the one-dimensional (ID) 



Bose gas in a harmonic trap of energy E = ^jrijtj, 
where single-particle energies ej — hwj In order to 
obtain the probability distribution we first consider par- 
tition functions. The canonical ensemble partition func- 
tion for a system with N^x excited bosons and tempera- 
ture T of ID harmonic oscillator is [lo| : 
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where nj are populations of jth state, Nex = J^JLi is 
a number of bosons and /? = l/fcsT. Kronecker delta 
function in ^ enforces the condition of N^x thermal 
bosons. In terms of these functions probability distri- 
bution of having N^x thermal bosons is P-^^ {N^x, P) — 

Zex{Nex,f3)/Z]\f with Zn = J^Nex {Nex,f3). Intro- 

ducing an integral representation of the delta function 
Sa,b = 57 /g e™'''~''''da: and calculating the integral we 
obtain: 
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Similar probability distribution can be calculated for the 
classical counterpart. In this case the canonical ensem- 
ble partition function must take into account continuous 
values of populations while energy levels are quantized 
as previously. For the classical counterparts sums are re- 
placed by integral and the canonical ensemble partition 
function is: 

ZexiNex, f3)= J^--- e-P'5(M,x - N,x) , 

(3) 

where amplitudes \aj\'^ G (0,oo), the energy £ = 

J2f=i'' \ot]V^i ^'^^ '7^ ex = X)j^T^ Notice, that 

we include a finite numbers of states up to Kmax only. 
Populations of the higher energy states are set to zero 
and not taken into account within the classical approx- 
imation. The integral representation of the Dirac delta 
function S{a — b) = ^ da:e*'^(°~''^ simplifies calcula- 
tions which give: 
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FIG. 1: (Color online) Probability distribution of having Nex 
bosons outside the ground state of a ID trap. Points repre- 
sent the exact quantum distribution while lines are its clas- 
sical counterparts with hio/SKmax = 0.9 - dot-dashed line, 
huPKmax = 1 - solid line and hujpKmax = 1.1 - dashed line. 
Total number of bosons is A'' = 1000 and at TujjP = 0.04. 
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FIG. 2: (Color online) Solution of the equation l|10p as a 
function of T/Tc for 3D harmonic trap and total number of 
bosons iV = 10^ (T), iV = 2 • 10'' (□), = 6 • 10^ (•). Dashed 
line represents a value of hujKmaxP given by the expression 



Note, that this formula leads to the equipartition of en- 
ergy in every single-particle mode, i.e. {nj)ej — ksT, 
where {rij) is the mean occupation of the state of energy 
Ej. This feature is characteristic for thermal state of a 
classical system. Probability distributions P^^ {N^x, P) 
and its classical counterpart V^^iNex, (3) are shown in 
FigdJ One can see that there exists an optimal value of 
the cut-off which in this case is equal to Kmax — l/huj(3. 
For this optimal value quantum distribution and its clas- 
sical counterpart are almost indistinguishable. Notice, 
that other values of Kmax lead to a shift of the entire 
distribution while its shape does not change significantly. 

It is interesting to analyze the system in 2 and 3 spatial 
dimensions. Now, the problem is more complex since we 
have to deal with a degeneracy. There exists excellent 
approximate formulas in the literature and we find the 
result of [1| the most useful for our purposes: 
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Zn {U/ri-l)\{N,x)\ 
with parameters: 
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Please note again that each term in the sum ([6]) is a 
mean occupation of the single-particle level j which in 
the classical version |[8]) implied equipartition (fT.j)ej — 
ksT. The simple expression ^ is valid for any trapping 
potential which controls functions TC and t] through the 
single particle energy spectrum ej of a trap. Within the 
classical fields formula |[5|) is still valid but parameters Ti 
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FIG. 3: (Color online) Scaling of the cut-off parameter 
hiuKmaxP with dimensionality D. Points are given by formula 
(HH), line is the best linear fit: hujKmaxP = 0.477 + 0.562D. 



and 77 are different: 
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This classical version of Eq.(l5]) is not obvious. We 
checked the classical counterpart of ^ by comparing 
with 1(4]) in ID, and with the saddle points method |ll| re- 
sults in 3D. For the systems with small number of bosons, 
as considered numerically here, the maxima of distribu- 
tions coincide but there are some deviations in its width. 
The agreement improves for larger systems. 

In order to compare the quantum and classical coun- 
terpart of the probability distribution ^ we choose the 
cut-off Kjnax in such a way that the maxima of both dis- 
tributions coincide. The maximum of probability distri- 
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bution ^ is equal to H and equation for optimal K„ 
leads: 



(10) 



In FigO the solution of the above equation is plotted 
for 3D geometry and different total number of bosons. 
For large systems the value of hujKmaxP practically does 
not depend on T. Some deviations from this constant 
value can be observed for very low relative temperatures 
T/Tc only (Tc is a critical temperature) and this region 
decreases with increasing N. The value of fuoKmaxP sat- 
urates at: 



I, iOTD = 1, 

\ aD)(D-l){D-iy., i0TD>2, 

(11) 

where C(^) is the Riemann Zeta function This con- 
dition can be used as a criterion for the optimal choice 
of the maximal energy of classical states hojKmax- The 
condition ifTTj) defines the universal value of the cut-off 
Kmax which for large systems, as considered in the classi- 
cal fields method, depends only on temperature and space 
dimensionality. In Fig|4]we compare 3D quantum proba- 
bility distribution given by ([5]) with its classical counter- 
part where the Kmax was chosen according to (fTO| . We 
see that not only the maxima coincide but also the entire 
distributions are very similar. Therefore, by the appro- 
priate choice of one parameter only, the cut-off energy 
of classical modes huiKmax , we can accurately mimic the 
true quantum statistics by the statistics of the classical 
fields satisfying the equipartition of energy. This optimal 
cut-off scales nearly linearly with dimensionality D, see 
FigH 

The same analysis can be applied for the uniform 
system with periodic boundary conditions, i.e. cubic 
box of length L. The single-particle energy spectrum 
is then ek = ?i^k^/2m, where k = 2Tr{nx,ny,nz)/L 
and m = 0,±1,±2,.... In this case, the cut-off 
parameter is determined by the maximal wave-vector 
T<-max = '2Tr{nmax,nmax,nmax)/L and the condition JlO]) 
defines its value. For large systems /3/i^K^^^/2m very 
weakly depends on temperature similarly as for the har- 
monic potential. Its value saturates at /?fi.^K^(j^/2m = 
7r(C(3/2)/4)2 ~ 1.34, 0.68, 0.29 in 3D, 2D, ID respec- 
tively what was found numerically. In Fig. [4] the quan- 
tum probability distribution and its classical counterpart 
are shown for the optimal cut-off given by Eq. I|10p . One 
can see that not only the maximum but also a shape of 
both distributions are very close. 

Having matched the probability distributions it is of 
interest to calculate the average population (n^f^^^) of 
the final mode of the classical distribution. It can be 
obtained by considering probabilities p(na > /) of a mode 
q being occupied by at least I bosons |l3| : 
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FIG. 4: (Color online) Probability distribution of having N^x 
bosons outside the ground state of a 3D trap (a) and 3D box 
(6) potentials. Points are exact distribution while solid lines 
represent its optimal classical counterparts. Total number of 
bosons is A'' = 1000. Notice, for 3D box potential probability 
distributions are wider as compared to 3D harmonic trap. 



under the condition of Nex thermal bosons. The proba- 
bility of having exactly I bosons in mode q is pin^^ — I) = 
p{nq > I) — p{nq > I Then the average population 

of q state, with the condition of Nex thermal bosons, is 
{n-q)Ncai = '-P(^q ~ 0- Summation over all Nex, 

with weight p{Nex\N) = Zex{Nex, P)/Zn of having Nex 
atoms within TV, gives the formula for the average popu- 
lation of the q mode: 
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(13) 



For the classical distribution, where populations change 
in the continuous way, it is convenient to use a probabil- 
ity density instead. In this case, the probability of hav- 
ing at least I bosons in mode q is pi — e~^'''^'^Zex{Nex — 
I) / Zex{Nex) while the probability density of having ex- 
actly I bosons in mode q is p(nq = I) = —dpi/Al. Then, 
the final expression for the average population of q mode 
is: 
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FIG. 5: (Color online) Average populations {nK„ax) of the 
final mode Kmax Eq. lflSll (•), and for the classical distribu- 
tion Ea. l|14p (■), for 3D harmonic trap (a) and 3D box (b) 
potentials. Lines represent expected results in the thermody- 
namical limit. Number of atoms A'^ = 1000. 



Eq. lfT4|) is a continuous version of Eq. (fT3|) and based on 
the Euler-Maclaurin summation formula the leading term 
of its difference is approximately equal to {nq)ci — {riq) — 
0.5. Therefore, the average occupation of the K„iax mode 
within the classical fields method is larger than the quan- 
tum one. In FigElthe average populations {riKmaJ) and 
{nKmar^)ci for 3D harmonic trap (a) and 3D box poten- 
tial (b) are shown as functions of the temperature T/Tc. 
Summation of populations over all excited modes must 
give a number of excited atoms N^.^. in both cases that 
is why in the classical case the populations have to be 
larger. The values of the average occupations are very 
close to expected ones in the thermodynamical limit. 

In summary, we considered a classical counterpart of 
the probability distribution of the ideal Bose gas within 
the canonical ensemble. We found that there exists the 
optimal cut-off which allows to match the probability dis- 
tribution of the condensate population by its classical 
fields counterparts. We found the scaling of the cut-off 
energies with temperature and dimensionality. This be- 
ing a very simple model it nevertheless sheds some light 
on the optimal choice of the cut-off for weakly interacting 
bosons and more significantly on the very foundations of 
the method itself. We understood now much better that 
finite number of the classical fields can adequately repre- 
sent the statistical properties of the full quantum many 
body system of bosons. 

Acknowledgment The authors are grateful for use- 
ful discussion with M. Brewczyk. The authors acknowl- 
edge support by the Polish Government research funds 
for 2006-2009. 



[1] M. Planck, Verhandlungen der Deutschen Physikalischen 

Gesellschaft 2, 202-204 (1900). 
[2] A. Einstein, Sitzber. Kgl. Preuss. Akad. Wiss. 1, 3 (1925); 

A. Einstein, Sitzber. Kgl. Preuss. Akad. Wiss. 22, 261 

(1924). 

[3] S.N. Bose, Z. Phys. 26, 178 (1924). 

[4] M.H. Anderson, J.H. Ensher, M.R. Mattews, C.E. Wie- 
man, E.A. Cornell, Science 269, 198 (1995); K. B. Davis, 
M.-O. Mewes, M.R. Andrews, N.J. van Druten, D.S. Dur- 
fee, D.M. Kurn, W. Ketterle, Phys. Rev. Lett. 75, 3969 
(1995). 

[5] N. N. Bogoliubov, J. Phys. U. S. S. R. 11, 23 (1947). 

[6] Y. Kagan, B. V. Svistunov, Phys. Rev. Lett. 79, 3331 
(1997); M. J. Davis, S. A. Morgan and K. Burnett Phys. 
Rev. Lett. 87, 160402 (2001); A. Sinatra, C. Lobo, Y. 
Castin, Phys. Rev. Lett. 87, 210404 (2001); K. Coral, 
M. Gajda, K. Rzfizewski, Opt. Express 8, 92 (2001); 
A. Sinatra, C. Lobo, Y. Castin, J. Phys. B: At. Mol. 
Opt. Phys. 35, 3599-3631 (2002); M. J. Davis and S. A. 
Morgan, Phys. Rev. A 68, 053615 (2003); M. Brewczyk, 
P. Borowski, M. Cajda, K. Rzqzewski, J. Phys. B 37, 
2725 (2004); M. J. Davis and P. B. Blakie, J. Phys. A: 
Math. Gen. 38, 10259 (2005); M. Brewczyk, M. Gajda, 



K. Rzcjiewski, J. Phys. B 40, R1-R37 (2007). 
[7] L. Zawitkowski, M. Brewczyk, M. Gajda, K. Rz^zewski, 

Phys. Rev. A 70, 033614 (2004); E. Witkowska, M. 

Gajda, J. Mostowski, J. Phys. B: At. Mol. Opt. Phys., 

40, 1-13 (2007). 
[8] We choose a zero of the energy scale at the single particle 

ground state. 

[9] V. V. Kocharovsky, VI. V. Kocharovsky, and Marian 
O. Scully, Phys. Rev. A 61, 053606 (2000); V. V. 
Kocharovsky, VI. V. Kocharovsky, and Marian O. Scully, 
Phys. Rev. Lett. 84, 2306 (2000). 

[10] C. Weiss and M.Wilkens, OPTICS EXPRESS 1 No. 10, 
272 (1997); S. Grossmann, M. Holthaus, OPTICS EX- 
PRESS 1 No. 10, 262 (1997); M. Holthaus and E. Kali- 
nowski, Ann. Phys. 270, 198-230, (1998). 

[11] M. Holthaus and E. Kalinowski, Ann. Phys. 276, 321- 
360, (1999). 

[12] Formula ifTTjl is easily derived with the help of the usual 
replacement of the sum by the integral in ^ and ^ and 
retaining the leading term in the degeneracy factor. 

[13] Z.Idziaszek, L. Zawitkowski, M. Gajda, and K. 
Rz^zewski, to be published. 



